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ABSTRACT 

We present a timing analysis of the 2002 outburst of the accreting millisecond pulsar 
SAX J1808.4— 3658. A study of the phase delays of the entire pulse profile shows a behavior that 
is surprising and difficult to interpret: superposed to a general trend, a big jump by about 0.2 in 
phase is visible, starting at day 14 after the beginning of the outburst. An analysis of the pulse profile 
indicates the presence of a significant first harmonic. Studying the fundamental and the first harmonic 
separately, we find that the phase delays of the first harmonic are more regular, with no sign of the 
jump observed in the fundamental. The fitting of the phase delays of the first harmonic with a model 
which takes into account the observed exponential decay of the X-ray flux (and therefore of the mass 
accretion rate onto the neutron star) gives important information on the torque acting on the neutron 
star during the outburst. We find that the source shows spin-up in the first part of the outburst, 
while a spin-down dominates at the end. From these results we derive an estimate of the neutron star 
magnetic field strength. 



Subject headings: stars: neutron — stars: magnetic fields 
SAX J1808.4-3658 — X-ray: binaries 



pulsars: general — pulsars: individual: 



1. INTRODUCTION 

In 1998, the idea that neutron stars in Low Mass 
X-ray Binaries (hereafter LMXBs) are spinning at mil- 
lisecond periods was spectacularly demonstrated by the 
discovery of coherent X-ray pulsations at ~ 2.5 ms in 
SAX J1808.4— 3658, a transient X-ray source with an or- 
bital period of 2 hr (Wijnands & van der Klis 1998). For 
almost four years SAX J1808.4— 3658 has been consid- 
ered as a unique object in which some peculiarity of the 
system allowed the detection of the neutron star spin. 
However, in the last few years six other accreting mil- 
lisecond pulsars have been discovered (see Wijnands 2005 
for a review): all of them are transient, with spin periods 
in the range between 1.7 and 5.4 msec. 

SAX J1808.4— 3658 was soon recognized as a weakly 
magnetized (< 10 10 G), rapidly rotating (401 Hz) ac- 
creting neutron star (Wijnands & van der Klis 1998; 
Chakrabarty & Morgan 1998). It is usually found in qui- 
escence, where it shows a luminosity below 10 32 ergs/s, 
and exhibits quasi-periodic outbursts, roughly every two 
years, with peak luminosity around 10 36 ergs/s. It shows 
all the phenomenology of a typical LMXB (i.e. type-I 
X-ray bursts and burst oscillations, kHz QPOs, etc.), to- 
gether with coherent pulsations, which makes this system 
unique. Although other X-ray millisecond pulsars have 
been discovered recently, SAX J1808.4— 3658, the only 
accreting millisecond pulsar for which more than one X- 
ray outburst has been observed, still remains the richest 
laboratory for the study of old accreting neutron stars in 
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binary systems. 

Still, no detailed timing studies have been published 
for this source so far. During the 2002 outburst of 
SAX J1808.4— 3658, the persistent accretion-powered X- 
ray pulsations were detected with a fractional r.m.s. am- 
plitude of 3-10%. The pulsar spin frequency derived from 
these data was approximately 400.97521 Hz at the start 
of the outburst, and a mean spin-down rate of about 
2 x 10 -13 Hz s" 1 has been claimed (Chakrabarty et al. 
2003, see also Morgan et al. 2003). 

In this paper we report the results of the detailed tim- 
ing analysis of the 2002 outburst of SAX J1808.4-3658; 
we find that the source shows a puzzling behavior of the 
phase delays, which can be interpreted analysing sepa- 
rately the fundamental and first harmonic of the pulse 
profile. 

2. TIMING ANALYSIS AND RESULTS 

SAX J1808.4-3658 was in outburst in 2002 between 
October 15 and November 26, and was extensively ob- 
served by Rossi X-ray Timing Explorer (RXTE). We 
analyse here all the available observations of the RXTE 
data archive taken during this outburst. We mainly use 
data from the RXTE Proportional Counter Array (PCA, 
Jahoda et al. 1996), which consists of five identical gas- 
filled proportional counter units (PCUs), with a total ef- 
fective area of ~ 6000 cm 2 , sensitive in the energy range 
between 2 and 60 keV. We used data collected in generic 
Events mode, with a time resolution of 125 /zs and 64 
energy channels. These files were processed and ana- 
lyzed using the FTOOLS v. 5. 3.1. In order to eliminate 
the Doppler effects caused by the Earth and satellite mo- 
tion, the arrival times of all the events were converted to 
barycentric dynamical times at the Solar system barycen- 
ter. The position adopted for the source was that of the 
proposed radio counterpart (0".4 uncertainty, which is 
compatible with that of the optical counterpart, Rupen 
et al. 2002; Giles, Hill & Greenhill 1999). The PCA 
light curve of the 2002 outburst of SAX J1808.4-3658 
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is shown in Figure 1 (upper panel). For the spectral 
analysis we also used data from the High-Energy X-ray 
Timing Experiment (HEXTE, - 20 - 200 keV energy 
range, Rothschild et al. 1998). 

To analyse these data we used the procedure that is 
extensively described in Burderi et al. (2006). We cor- 
rected the arrival times of all the events for the delays 
caused by the motion of the neutron star in the binary 
system, using the most updated orbital parameters (Pa- 
pitto et al. 2005). To compute phases of good statistical 
significance we divided the whole observation in intervals 
of 2 h, approximately the orbital period of the system: 
in this way we average over periodic phase residuals, if 
any, caused by errors in our assumed orbital parame- 
ters (see e.g. Galloway et al. 2002; Burderi et al. 2006). 
We epoch folded each interval of data in which the pul- 
sation was significantly detected at the spin period of 
2.493919760956 ms with respect to the same reference 
epoch, Xb = 52562 MJD, corresponding to the beginning 
of the outburst. The fractional part of the phase was 
obtained by fitting each pulse profile with two sinusoids 
of fixed periods (corresponding to the fundamental, with 
period equal to the spin period, and the first harmonic, 
with period equal to half the spin period, respectively). 
In this way we obtain the phase delays shown in Figure 2 
for the fundamental (right panel) and first harmonic (left 
panel), respectively. 

The phase delays of the fundamental show a puzzling 
behavior, with a clear jump by about 0.2 in phase which 
occurs at day 14 after the beginning of the outburst. In- 
terestingly, in the outburst light curve, the 14th day cor- 
responds to a change of the steepness of the exponential 
flux decay (see Fig. 1, upper): while before this day the 
X-ray flux decreases exponentially with a characteristic 
time of about 10 days, after this day the characteris- 
tic time of the exponential decay becomes about 3 days 
(a similar behavior of the light curve was observed for 
the 1998 outburst of SAX J1808.4-3658, Gilfanov et al. 
1998). After day 17 from the beginning of the outburst 
the X-ray flux shows rapid oscillations of the count rate 
on timescales of hours to days. On the other hand, the 
phase delays of the first harmonic do not show any evi- 
dence of such a jump. Note that this is not an effect of 
the larger error bars of the phases of the harmonic, as it 
can be seen from Figure [T] (lower panel), where we plot 
the differences between the phases of the fundamental 
and the phases of the harmonic as a function of time. 
For a global shift of the pulse profile the phases of the 
fundamental and harmonic should shift in the same way 
and the phase differences should remain constant. On 
the other hand, a variation of the differences indicates a 
shift of the phase of the fundamental with respect to the 
harmonic caused by a change of the shape of the pulse 
profile. In the plot of the phase differences it is possible 
to see that the big jump at day 14 is still visible in the 
differences, and is much larger than the associated error 
bars, meaning that it is present in the fundamental and 
not in the harmonic. We conclude that in this case the 
behavior of the harmonic is genuinely simpler than the 
behavior of the fundamental, and therefore we tried to 
fit it to an appropriate model. 

In order to derive the differential correction to the spin 
frequency, A^o, and the spin frequency derivative, v, 
we have to derive a functional form for the time de- 



pendence of the phase delays. For a constant varia- 
tion of the spin frequency (v = const), the phase de- 
lays (which are defined as the double integration over 
time of v) will be a parabolic function of time. However, 
since for SAX J1808.4-3658 the X-ray flux is observed 
to vary with time, we expect a more complex dependence 
of time. To derive this expression, we use the following 
simple assumptions: i) The bolomctric X-ray luminosity 
L is a good tracer of the mass accretion rate, M, via 
the relation L = (GM/R)M, where G, M, and R are 
the Gravitational constant and the neutron star mass 
and radius, respectively. Therefore, we assume for the 
mass accretion rate the following dependence on time: 
M(t) — Moexp[(t — T )/t], where r, the characteristic 
time of the flux exponential decay, that is 9.27 days from 
a fit of the X-ray light curve, ii) The accreted matter 
transfers to the neutron star its specific angular momen- 
tum at the corotation radius (see e.g. Rappaport et al. 
2004). In general, the matter should transfer to the neu- 
tron star its specific angular momentum at the accre- 
tion radius (that is the magnetospheric radius, the ra- 
dius at which the Keplerian accretion disk is truncated 
by the neutron star magnetic field). Since for accreting 
millisecond pulsars, the magnetospheric radius must be 
very close to the corotation radius (about 3 neutron star 
radii in the case of SAX J1808.4— 3658), our assump- 
tion can be considered a good approximation, and, in 
any case, gives a lower limit to the mass accretion rate 
needed to obtain the observed torque. Therefore, the rate 
of angular momentum transferred to the neutron star is 
L = 2tx1v = M(GMR C o) 1/2 , where I is the neutron 
star moment of inertia, v its spin frequency, and Rco 
the corotation radius. We do not consider here any form 
of threading of the accretion disk by the magnetic field 
of the neutron star (see e.g. Rappaport et al. 2004). 

From these assumptions, with some algebra, we derive 
the following expression for the phase delays caused by 
accretion torque due to an M decreasing exponentially 
with time: 

cf>(t) =(/) - B(t - T ) - Cexp[-(i - T )/r] (1) 

where (f>o is a constant, all the times are expressed in 
days, C = 1.067 x lO" 4 /^ 1 P% 3 m 2 / 3 T 2 m , where J 4S is 
the moment of inertia in units of 10 45 g cm 2 , P_3 is the 
spin period in msec, m is the neutron star mass in Solar 
masses, and rho is the accretion rate at t = Tq in units 
of 10~ 10 M yr" 1 , and, finally, B = Au + C/t, where 
Afo is the differential correction to the spin frequency 
adopted for the folding of the light curves. Note that 
we are neglecting here the change to a steeper fall off of 
the light curve which occurs at day 14. However, as it 
will become clear in the following, after day 14 the mass 
accretion rate is already so low that any spin-up torque 
becomes negligible (with respect to the behavior of the 
phase delays which instead tend to flattens, see below). 
Therefore, including a reduction of r in the model after 
day 14 (which further reduces the already negligible spin- 
up), does not change the results. 

We tried to fit the phase delays of the first harmonic 
with eq. ([T]), but we obtained a poor description of the 
data, corresponding to a \ 2 /dof = 613. 1/199. 5 Indeed, 

5 Note that fitting these data to simple linear or parabolic 
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using this model, we can obtain a good fit of the first 
14 days of the outburst, but, with respect to this fit, we 
observe a flattening of the phase delays after day 14. To 
describe this flattening we added to the model described 
by eq. ([1]) a quadratic term corresponding to a constant 
spin down. In this way we obtained a significant im- 
provement of the fit (x /dof = 485.0/198; this x 2 is still 
large, but now we do not see any systematic trend in the 
residuals with respect to the best fit model). Although 
in this way we have obtained a great improvement of 
the quality of the fit, the obtained x 2 is still unaccept- 
able, because of localised residuals around day 4 and day 
30 — 31 from the beginning of the outburst. These resid- 
uals cannot be taken into account modifying the fitting 
function, but the uncertainties derived from a fit that is 
not acceptable in the statistical sense may be underesti- 
mated. To obtain more conservative uncertainties for our 
best fit parameters, we therefore increased by a factor 1.5 
the errors of all our phase points, obtaining in this way 
a x 2 /dof — 215/198 (that is now very close to 1), and 
we have re-calculated all the uncertainties on our best fit 
parameters. The best fit parameters and the correspond- 
ing uncertainties (at 90% c.l.) are shown in Table 1; the 
best fit function is plotted on top of the data in Figure 2. 
For a comparison, we have plotted the best fit function 
of the first harmonic on top of the fundamental. 

The last issue we have to discuss regards the sys- 
tematic uncertainties induced by the uncertainty on the 
source position. As discussed in Burderi et al. (2006), 
these systematics give both a systematic uncertainty 
on the linear term (that is the spin frequency correc- 
tion) and on the quadratic term (that is the spin fre- 
quency derivative). For an error circle of 0".4 around 
the position of SAX J1808.4-3658, we find that the sys- 
tematic uncertainty on this two terms are, respectively, 
Ai/ sys = 7.9 x 10" 8 Hz and Ai> sys = 1.8 x 10" 14 Hz 
s — 1 . This increases the total uncertainty on the spin fre- 
quency to 1.1 x 10~ 7 Hz, and the total uncertainty on 
the spin-up and spin-down terms, which become: vq — 
(4.40±0.83) x 10~ 13 Hz s" 1 and ;> sd = -(7.6±2.3) x 10~ 14 
Hz s _1 , respectively. Even considering these systematic 
uncertainties, both the spin-up and spin-down terms are 
still significant at more than 4cr (note that all the re- 
ported uncertainty are at 90% c.l.). 

Finally, in order to compare the derived spin-frequency 
derivative with the mass accretion rate Mq as inferred 
from the bolometric X-ray flux of the source, we per- 
formed a spectral analysis of the PCA and HEXTE spec- 
tra at the beginning of the outburst (i.e. at Oct 15 and 
16, 2002). Keeping the absorption to the source fixed at 
Nh = 3.0 x 10 21 cm~ 2 (the Galactic X-ray absorption 
in the direction of SAX J1808.4— 3658 at a distance of 
3.5 kpc, Galloway & Cumming 2006), a very good fit 
to the 2.5 - 200 keV spectrum of SAX J1808.4-3658 is 
given by a disk blackbody plus a cutoff power law, and a 
Gaussian emission line with centroid fixed at 6.4 keV (the 
fluorescence K-shell iron line). We extrapolated the ob- 
served X-ray flux in a broad band energy range in order 
to evaluate the bolometric accretion luminosity, which re- 
sulted to be ~ 1 x 10 37 ergs/s on 2002 Oct 15 assuming 
a distance of 3.5 kpc and correcting for the interstellar 

functions gives similar (slightly worse) results, corresponding to 
X 2 /dof = 623.4/200 and X 2 /dof = 622.4/199, respectively. 



absorption. 

3. DISCUSSION AND CONCLUSIONS 

The results described in the previous section demon- 
strate that the timing analysis for an accreting pulsar is 
more complex than expected. The phase delays of the 
fundamental clearly show a jump by 0.2 in phase which 
may be due to instabilities induced by the accretion of 
matter onto a weakly magnetized star. Indeed, if the 
magnetic field is not strong enough to completely dom- 
inate the motion of matter to the polar caps, we can 
expect that variations in the accretion flow may cause 
small movements of the footpoint on the neutron star 
which can give rise to the shift in phase that we ob- 
serve. On the other hand the phase delays of the first 
harmonic show a more regular behavior and can be fit- 
ted to a model derived from the theory of the accretion 
torque, as it will become clear in a while. Therefore, 
the first harmonic might have some fundamental physi- 
cal meaning, probably related to the fact that this may 
represent the accretion onto both the polar caps. 

As mentioned above, the results of the fitting of the 
phase delays of the harmonic are in agreement with 
torque onto the neutron star predicted by the accretion 
theory. A simple second order polynomial or a model 
taking into account the exponential decrease with time 
of M do not give a good fit of the phase delays in the 
whole time range; these models give a good fit if we con- 
sider only the first 15 days of the outburst. We obtain 
a good fit of the phase delays during the whole outburst 
using the model described by eq. Q] and adding to this 
model a quadratic term, which describes the flattening 
of the phases at the end of the outburst. This means 
that the neutron star is spinning-up at the beginning of 
the outburst, as expected when the mass accretion rate 
is relatively high, but spins-down at the end of the out- 
burst. 

This gives very important information on the torque 
acting onto the neutron star in this system, implying 
that, when the mass accretion rate has significantly de- 
creased, the torque onto the neutron star changes sign. 
This can be caused, for instance, by a threading of the 
accretion disk which becomes important at low accretion 
rate (see eq. (23) in Rappaport et al. 2004). In this case, 
we can evaluate the magnetic field of SAX J1808.4— 3658 
from our measured value of the spin-down using the re- 
lation: /i 2 /(9r^ Q ) = 27r/i> s d. The magnetic field found 
in this way is B ~ (3.5 ± 0.5) x 10 s Gauss, perfectly in 
agreement with previous constraints (see e.g. Psaltis & 
Chakarabarty 1999; Di Salvo & Burderi 2003). 

On the other hand, the spin-up observed in the first 
part of the outburst will give information on the mass 
accretion rate at the reference time t = To, that is 
Mo — 1-8 x 10~ 9 M yr _1 . This gives a bolometric 
X-ray luminosity at the beginning of the outburst of 
~ 2 x 10 37 ergs/s. This is about a factor of 2 higher 
than the bolometric X-ray luminosity inferred by the 
SAX J1808.4-3658 RXTE spectra (about 1 x 10 37 ergs/s 
assuming a distance of 3.5 kpc). On the other hand, as- 
suming that the mass accretion rate, Mo, inferred from 
our timing analysis is correct, from a comparison with 
the observed X-ray flux, we can conclude that the source 
should be at a distance of ~ 5 kpc. 
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Table 1. Best fit parameters of the phase delays of the first har- 
monic. 



Parameter 


Value 


Parameter 


Value 


B 


(-0.8 ± 2.6) X 10~ 3 


Av 


-(3.61 ±0.72) x 10" 7 Hz 


C 


0.282 ±0.054 




(4.40 ± 0.81) x 10" 13 Hz s" 1 


D 


(2.85 ±0.59) x 10" 4 




-(7.6 ± 1.5) x 10~ 14 Hz s- 1 


M 


(1.81 ± 0.32) X 10" 9 M yr - 1 


vo 


400.975209690 ± 7.2 X 10" 8 Hz 



Note. — The best fit function is given by eq. JTJ + D(t — To) 2 , where D = (1/2)!> s( j. The reference time To is the beginning of the 
observation, that is 52562 MJD. Quoted errors are at 90% c.l. and do not include the systematic uncertainties discussed in the text. 
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Time (days since 52562.0 MJD) 

Fig. 1.— Upper: PCA count rate of SAX J1808. 4-3658 during the 2002 outburst. The start time, T , is 52562 MJD, corresponding to 
October 15 2002. Lower: Phase differences, fundamental minus harmonic, vs. time. The jump of the phase of the fundamental at day 14 
is clearly visible in the differences. 



6 




10 20 30 40 10 20 30 

Time (days since 52562.0 MJD) Time (days since 52562.0 MJD) 



Fig. 2. — Left: Phase vs. time for the fundamental of the pulse frequency of SAX J1808.4— 3658. Right: Phase vs. time for the first 
harmonic of the pulse frequency of SAX J1808.4— 3658. On top of the data, the best fit function (the sum of eq. |(T} and a quadratic term) 
is plotted as a solid line. 



